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ABSTRACT 


Progress  in  research  into  fundamental  mechanisms  of  combustion  in  high 
speed  air  flows  is  described.  Three  major  areas  are  covered:  the  application 
of  a  modular  combustor  model  for  the  liquid  fueled  sudden  expansion  combustor 
configuration  to  the  analysis  of  flame  stability  and  the  interpretation  of 
combustor  test  data;  the  extension  of  the  modular  model  to  the  analysis  of 
boron  slurry  fueled  ramjet  combustors;  and  the  development  of  analytical  models 
of  the  ducted  rocket  combustor.  The  modular  approach  is  demonstrated  to  be 
a  useful  tool  for  the  analysis  and  interpretation  of  sudden-expansion  combustor 
flame  stabilization  phenomena  and  to  be  of  considerable  utility  in  the  analysis 
of  combustor  test  data.  Models  for  the  combustion  of  boron  particles  are 
reviewed  and  put  into  a  form  suitable  for  incorporation  in  a  modular  model  of 
a  boron  slurry  fueled  combustor.  Particle  tracking  methods  suitable  for  use 
in  the  context  of  a  modular  model  formulation  have  been  defined  and  developed 
for  the  required  application.  The  approach  required  to  provide  an  analytical 
model  of  the  three-dimensional  recirculating  flow  in  a  ducted  rocket  combustor 
has  been  defined  and  the  development  of  an  aerodynamic  model  has  been  initiated. 


1.  INTRODUCTION 


Design  trends  for  current  and  future  airbreathing  propulsion  systems  are 
trending  toward  lower  weight  and  smaller  volume  designs.  These  requirements, 
and  increased  demands  for  higher  performance  from  airbreathing  propulsion 
systems  in  general  and  ramjet  devices  in  particular,  have  resulted  in  the  need 
to  upgrade  existing  technologies  to  meet  the  emergence  of  more  stringent 
design  requirements.  Achieving  reduced  volume  and  weight  requires  the  use  of 
shorter  combustors  and  high  energy,  high  density  fuels.  Systems  constraints 
have  dictated  the  use  of  sudden-expansion  (dump)  combustors  capable  of 
operating  effectively  at  high  combustion  intensities  over  wide  ranges  of 
conditions.  Problems  of  flame  stabilization,  flame  propagation  and  spray 
combustion  have  been  encountered.  To  aid  in  the  development  of  solutions  to 
these  problems,  analytical  combustor  models  are  a  requirement:  the  ability 
to  compute  in  some  detail  combustion  chamber  flowfields  is  necessary  in  order 
to  understand  the  phenomena  that  occur  in  existing  combustors  and  to  predict 
the  performance  of  new  combustor  concepts.  The  insight  gained  through  the  use 
of  analytical  combustor  models  can  be  of  substantial  value  in  the  planning  of 
a  combustor  test  program  and  in  the  interpretation  of  combustor  and  combustor 
component  test  data. 

In  response  to  the  needs  just  outlined,  a  detailed  model  of  the  sudden- 
expansion  liquid  fueled  ramjet  combustor  has  been  developed,  as  part  of  earlier 
AFOSR-sponsored  work  (Ref.  1).  This  model  makes  use  of  the  modular  concept, 
in  which  the  combustor  flowfield,  represented  schematically  in  Fig.  1,  is 
broken  down  into  three  major  components:  a  directed  flow,  which  is  treated 
as  parabolic,  a  recirculation  zone,  assumed  to  be  represented  by  well- 
stirred  reactor(s),  and  a  turbulent  shear  layer  along  the  dividing  streamline 
which  separates  the  other  two  regions.  The  shear  layer  serves  as  the  coupling 
region  between  the  other  two  model  components;  fluxes  of  species  and  energy 
across  this  shear  layer  form  the  boundary  conditions  on  the  two  computational 
regions.  Finite-rate  chemistry,  based  on  the  quasiglobal  model  (Ref.  2)  is 
included  in  the  formulations  for  both  the  directed  flow  and  well-stirred 
reactor  regions ,  although  for  the  modular  model  calculations  described  in  this 
paper,  the  recirculation  region  well -stirred  reactor  formulation  has  been 
restricted  to  a  global  finite-rate  chemistry  model.  The  directed  flow  is 
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Figure  1.  Schematic  of  Sudden-Expansion  (Dump)  Burner. 
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assumed  to  be  fully  turbulent,  with  the  turbulent  viscosity  defined  by  a  two- 
equation  turbulence  model  (Ref.  3).  A  key  feature  of  this  approach  is  the 
provision  for  the  shear  layer  coupling  region  in  the  model.  Through  the  use 
of  this  element  of  the  model,  the  division  of  the  mass  flux  between  the 
directed  flow  and  the  recirculation  region  is  computed  iteratively  rather  than 
specified  empirically.  Furthermore,  the  directed  flow  region  is  computed  in 
detail  as  a  two-dimensional  parabolic  flowfield,  rather  than  through  a  one¬ 
dimensional  approximation,  allowing  the  use  of  detailed  computations  of  the 
mixing  and  chemical  reactions  in  this  region  of  the  combustor. 

A  variety  of  results  obtained  in  the  process  of  validating  the  modular 
model  approach  for  the  prediction  of  ramjet  combustor  performance  and  flow- 
field  characteristics  have  been  described  in  Refs.  1  and  4.  During  the 
current  program,  attention  has  been  focused  on  further  development  of  certain 
of  the  fundamental  elements  of  the  model,  with  particular  emphasis  on  the 
sti rred-reactor  representation  of  flame-stabilizing  recirculation  regions. 
Other  aspects  of  the  current  work  include  the  development  of  boron  combustion 
models  and  the  definition  of  modeling  approaches  for  ducted  rocket  combustor 
configurations. 

Both  gas  generator  ramjets  (ducted  rockets)  and  slurry-fueled  ramjets 
are  examples  of  propulsion  devices  which  involve  the  turbulent  mixing  and 
combustion  of  multi-phase  flows.  In  both  of  these  applications  boron-loaded 
fuels  are  of  considerable  interest  because  of  their  high  energy  content,  and 
for  both  of  these  systems  it  is  desirable  to  develop  models  which  can  aid 
in  obtaining  a  fundamental  understanding  of  the  phenomena  involved  in  two- 
phase  mixing  and  combustion  and  which  can  be  used  to  parametrically  assess 
combustor  performance.  This  need  is  being  addressed  in  two  of  the  tasks  that 
are  part  of  the  current  program.  For  boron  slurry-fueled  ramjets,  the 
approach  being  taken  involves  a  direct  extension  of  the  modular  modeling 
approach  developed  for  liquid-fueled  ramjets  and  described  in  Refs.  1  and  4, 
while  the  application  of  modular  modeling  to  the  more  complex  geometry  associ¬ 
ated  with  the  ducted  rocket  configuration  is  also  being  investigated.  In 
the  latter  work,  the  approach  being  followed  involves  the  development  of  a 
three-dimensional  elliptic  Navier-Stokes  solution  technique  to  be  used  to 
delineate  the  flowfield  regions  that  are  to  be  described  by  separate  elements 
of  the  modular  approach.  This  solution  procedure  is  being  developed  to 
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obtain  the  cold-flow  aerodynamics  of  the  combustor  configuration ;  the  detailed 
modeling  of  this  configuration  incorporating  the  effects  of  multiphase  flow 
and  chemical  heat  release  is  to  be  carried  out  using  the  modular  approach. 

The  three  subsequent  sections  of  this  report  reflect  the  three  main 
areas  of  investigation  outlined  above.  Thus,  in  Section  2,  further  examination 
of  the  use  of  the  modular  model  of  a  liquid-fueled  sudden  expansion  combustor 
is  described,  with  emphasis  on  the  prediction  of  flame  stabilization  phenomena. 
A  detailed  model  of  the  combustion  dynamics  inherent  in  a  boron  slurry-fueled 
ramjet  combustor  requires  the  development  of  models  for  the  combustion  of  boron 
particles  and  of  techniques  for  treating  the  details  of  the  combustion  of 
boron  particle  clouds,  and  progress  in  this  area  is  discussed  in  Section  3. 
Then,  in  Section  4,  the  technique  being  utilized  to  develop  a  gas-generator 
fueled  ramjet  (ducted  rocket)  aerodynamic  model  is  described,  and  the  overall 
progress  made  during  this  program  is  summarized  in  Section  5.  In  each  of  these 
areas  close  coordination  has  been  maintained  with  government  laboratory  experi¬ 
mental  work,  including  liquid-fueled  ramjet  and  ducted  rocket  investigations 
underway  at  the  Aeropropulsion  Laboratory  of  the  Air  Force  Wright  Aeronautical 
Laboratory  (AFWAL)  and  boron-slurry  ramjet  and  ducted  rocket  work  underway 
at  the  Naval  Weapons  Center  (NWC)  under  joint  AFOSR/NWC  sponsorship. 


2.  THE  APPLICATION  OF  MODULAR  MODELING  TO  RAMJET 


PERFORMANCE  INVESTIGATIONS* 


The  well -stirred  reactor,  shown  schematical ly  in  Fig.  2a,  is  a  laboratory 
device  in  which  very  high  mixing  rates  are  achieved.  In  general,  laboratory 
stirred  reactors  are  designed  to  ensure  that  the  mixture  within  the  reactor 
is  spatially  uniform,  so  that  unreacted  feed  material  is  continuously  and 
uniformly  mixed  with  combustion  products,  reacting  for  a  time  defined  by  the 
average  residence  time  of  the  reactor  before  exiting.  For  the  limit  of  perfect 
stirring,  this  average  residence  time  is  given  simply  by  the  ratio  of  the 
stirred  reactor  volume  to  the  feed  mass  flux,  V/m;  the  product  mass  flux  is 
of  course  equal  to  the  feed  mass  flux.  In  this  limit  the  governing  equations 
for  the  stirred  reactor  state  reduce  to  algebraic  relations  (Ref.  4),  which 
allows  the  use  of  rapid  and  efficient  solution  procedures,  incorporating  either 
global  finite-rate  kinetics  or  the  more  detailed  quasiglobal  kinetics  formulation. 

In  the  lower  half  of  Fig.  2,  a  typical  recirculation  region  in  a  sudden- 
expansion  combustor  is  sketched.  Because  of  the  high  turbulence  intensity  and 
large  mixing  rates  generated  in  the  flow  reversal  within  the  recirculation 
the  limiting  behavior  of  this  region  can  be  thought  of  as  well-mixed;  that  is, 
the  state  of  the  recirculation  region  is  defined  by  the  chemical  kinetics  of  the 
reactions  occurring  and  not  by  the  mixing  rate.  In  this  sense  the  volume  within 
the  recirculation  region  sketched  in  Fig.  2b  is  similar  to  the  volume  within  the 
laboratory  stirred  reactor  of  Fig.  2a,  and  the  same  solution  technique  can  be 
used  to  obtain  the  thermochemical  state.  A  comparison  of  Figs.  2a  and  2b  also 
shows  the  major  difference  between  the  recirculation  region  and  the  laboratory 
stirred  reactor.  In  the  laboratory  stirred  reactor  discrete  reactant  and  prod¬ 
uct  streams  can  be  identified,  but  in  the  well-stirred  reactor  model  of  a 
sudden-expansion  recirculation  zone,  separate  reactant  and  product  streams  cannot 
be  defined.  Instead,  reactants  enter  the  recirculation  region  and  products  leave 
the  recirculation  region  by  turbulent  diffusion  through  the  shear  layer,  inte¬ 
grated  over  the  surface  area  of  the  dividing  streamline,  and  the  product  stream 
is  defined  similarly.  These  fluxes  are  of  course  equal  and  opposite. 


*  The  material  incorporated  in  this  section  is  taken  from  "Interpretation  of 
Ramjet  Combustor  Test  Data,"  by  P.  T.  Harsha  and  R.  B.  Edelman,  AIAA  Paper 
81-1433,  to  be  presented  at  the  AIAA/SAE/ASME  17th  Joint  Propulsion  Conference, 
July  27-29,  1981. 
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(a)  Schematic  design  of  spherical 
well-stirred  reactor. 


(b)  Schematic  of  sudden-expansion 
recirculation  region. 


FIGURE  2.  Comparison  of  well-stirred  reactor 
with  recirculation  region  flowfield. 
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The  fact  that  in  the  well -stirred  reactor  model  of  a  sudden-expansion 
recirculation  region  the  feed  rates  are  defined  by  the  fluxes  of  reactants 
and  products  through  the  shear  layer  region  is  the  reason  that  this  region 
has  been  identified  as  a  key  element  of  the  modular  model  described  in  Refs. 

1-3  and  6.  Nevertheless,  this  shear  layer  is  modeled  fairly  crudely  as  a  region 
of  linear  gradients:  that  is,  the  gradient  in  a  quantity  <$>  at  the  dividing 
streamline  is  approximated  by 


where  the  subscripts  P,  R,  and  W  refer  to  conditions  in  the  directed  flow 
immediately  outside  of  the  recirculation  zone,  in  the  recirculation  zone  itself, 
and  along  the  dividing  streamline,  respectively.  Moreover,  the  width  of  the 
shear  layer  is  assumed  to  be  given  by  the  linear  expression 

£  =  ax  +  b  (2) 

in  which  a  and  b  are  constants  to  be  defined  through  comparison  of  model 
predictions  with  experimental  data. 

As  is  described  in  more  detail  in  Refs.  1,  4-6,  the  overall  flowfield 
computation  using  the  modular  model  proceeds  as  follows:  a  dividing  stream¬ 
line  shape  is  assumed,  and  the  shear  layer  width  expression  and  shear  stress 
distribution  along  the  dividing  streamline  are  specified,  the  latter  through 
use  of  a  specified  "skin  friction"  coefficient.  An  initial  state  for  the 
stirred  reactor  is  assumed,  which  defines  an  initial  guess  for  the  species 
distribution  and  temperature  within  the  recirculation  zone.  With  these 
available  as  boundary  conditions,  a  finite-difference  computation  of  the  mixing 
region  external  to  the  recirculation  is  carried  out  to  the  end  of  the  recircu¬ 
lation  region.  This  calculation  defines  the  species  mass  fraction  and 
temperature  gradients  along  the  dividing  streamline,  and  these  values  are  used 
as  outlined  above  to  obtain  the  stirred  reactor  feed  rates.  A  new  stirred 
reactor  computation  is  carried  out  using  these  feed  rates,  resulting  in  a  new 
specification  of  the  species  mass  fractions  and  temperature  within  the  recircu¬ 
lation  region,  and  the  parabolic  mixing  computation  is  repeated.  This  procedure 
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is  continued  until  changes  in  the  stirred  reactor  composition  and  temperature 
from  computation  to  computation  become  small,  typically  0.1%,  at  which  point 
the  coupling  iteration  has  converged.  The  parabolic  calculation  is  then 
carried  out  to  the  end  of  the  combustion  chamber,  completing  the  solution. 

The  entrainment  rate  into  the  recirculation  region  is  an  important 
parameter  in  the  overall  combustor  flowfield  computation,  and  as  the  previous 
paragraph  indicates,  this  parameter  is  not  specified,  but  is  calculated  as  a 
part  of  the  modular  model  solution  procedure.  Thus  it  is  of  considerable 
interest  to  compare  the  entrainment  rate  computed  by  the  model  with  experimental 
data.  A  variety  of  data  for  recirculation  region  entrainment  rates  was 
reviewed  by  Curran  (Ref.  7)  and  the  range  of  the  observed  correlations  is 
shown  in  Fig.  3:  these  correlations  are  for  both  cold  flow  and  hot  flow. 

Also  shown  in  Fig.  3  are  the  entrainment  mass  flux  ratios  predicted  using  the 
modular  model  in  several  different  calculations,  involving  different  area 
ratios  and  cold  flow  as  well  as  reacting  flow.  Other  results  from  many  of 
these  computations  are  described  in  Ref.  4.  The  range  of  correlations  shown 
did  not  include  data  from  the  investigations  for  which  the  calculations  shown 
in  Fig.  3  were  made;  nevertheless,  the  predictions  show  reasonably  good  agree¬ 
ment  with  the  observed  correlations.  Of  particular  significance  is  the  fact 
that  the  predicted  reduction  of  the  entrainment  rate  with  combustion  compared 
to  cold  flow  at  the  same  area  ratio  is  in  agreement  with  the  observations 
reported  by  Curran.  This  represents  an  important  result,  both  from  the 
standpoint  of  validation  of  the  modular  model  approach,  even  though  a  rela¬ 
tively  crude  shear  layer  model  has  been  used,  and  because  it  suggests  a  more 
systematic  study  that  can  be  pursued  using  the  analytical  model  to  examine 
a  wider  range  of  operating  conditions  than  those  on  which  the  correlations 
were  originally  based. 

One  of  the  most  important  determinants  of  liquid  fueled  ramjet  performance 
is  the  fuel  distribution  at  the  combustor  entrance.  For  the  sudden-expansion 
liquid  fueled  ramjet,  the  interaction  between  inlet  fuel  distribution  and 
recirculation  zone  state  is  a  further  important  consideration.  A  variety 
of  observations  have  shown  that  for  wall  injection  of  the  fuel,  the  recircula¬ 
tion  zone  equivalence  ratio  is  generally  more  fuel-rich  than  the  overall 
(global)  equivalence  ratio.  This  is  demonstrated  by  the  data  of  Schmotolocha 
and  Economos  (Ref.  8),  which  was  obtained  for  three  different  inlet  conditions 


8 


IGURE  3.  Comparison  of  measured  and 

nredicted  entrainment  rates. 


and  at  two  different  locations  of  the  fuel  injection  orifices  upstream  of  the 
dump  plane.  Indeed,  the  measurements  described  in  Ref.  8  indicate  that  for 
overall  equivalence  ratios  greater  than  0.2,  the  equivalence  ratio  in  the 
recirculation  region  is  always  greater  than  unity.  These  data  are  shown  in 
Fig.  4,  in  which  the  results  of  a  modular  model  computation  at  one  overall 
equivalence  ratio  are  also  shown.  As  was  the  case  for  the  entrainment  rate 
pedictions,  the  agreement  between  the  experimental  results  and  the  computation 
is  reasonably  good. 

In  the  experiments  described  in  Ref.  8  plain  orifice  fuel  injection  was 
used.  For  plain  orifices,  fuel  penetration  characteristics  are  reasonably 
well  documented,  and  an  empirical  penetration  correlation  was  used  to  establish 
initial  conditions  for  the  computations  used  to  obtain  the  result  shown  in 
Fig.  4.  However,  plain  orifices  are  not  widely  used  in  ramjet  design  because 
of  their  poor  turndown  characteristics  which  result  in  the  need  for  several 
sets  of  orifices  to  obtain  a  wide  ramjet  operating  range.  The  most  commonly 
used  fuel  injector  in  liquid  fuel  ramjet  applications  is  the  poppet  injector, 
which  has  a  much  broader  turndown  range  than  does  a  plain  orifice.  However, 
penetration  characteristics  from  poppet  injectors  under  ramjet  operating 
conditions  are  not  well  documented.  Indeed,  there  exists  one  investigation 
of  poppet  (and  other  nozzle)  penetration  characteristics  under  ramjet  operating 
conditions  in  which  it  was  observed  that  the  penetration  from  poppets  is 
essentially  nil.  Ref.  9.  While  other  results  indicate  that  penetration 
characteristics  similar  to  orifice  injection  can  be  obtained  using  poppets. 

Ref.  10,  it  remains  clear  that  the  penetration  characteristics  of  poppets 
are  not  as  well  known  as  for  plain  orifices. 

At  one  operating  condition,  which  involves  a  relatively  high  inlet 
temperature,  poor  ramjet  combustor  performance  was  observed  in  a  recent 
development  program:  at  this  one  condition  performance  was  considerably 
poorer  than  at  other  test  conditions.  Review  of  the  available  results  suggested 
that  fuel  distribution  effects  were  a  possible  cause,  and  in  particular  that 
at  the  operating  condition  in  question  the  fuel  penetration  from  the  wall- 
mounted  fuel  injectors  may  have  been  negligible.  To  investigate  this  possibility, 
and  to  study  the  effects  of  fuel  distribution  on  combustor  performance,  a  series 
of  parametric  modular  model  computations  was  carried  out. 
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In  the  modular  approach,  fuel  injection  is  modeled  using  an  empirical 
correlation  for  penetration  and  breakup  time  correlation.  The  latter  correla¬ 
tion  defines  the  downstream  position  at  which  the  penetration  is  computed; 
from  this  point  the  fuel  spray  is  assumed  to  vaporize  at  a  rate  given  by  a 
bulk  spray  vaporization  correlation.  After  injection,  the  fuel  is  assumed 
to  reside  in  an  annular  region  at  the  position  specified  by  the  penetration 
correlation,  and  spreading  of  this  annulus  is  computed  using  a  turbulent  mixing 
model.  For  the  conditions  of  interest  in  the  parametric  study,  vaporization 
occurs  rapidly,  so  that  the  parameter  to  be  varied  was  the  fuel  penetration 
height. 

Fig.  5  shows  the  predicted  recirculation  zone  equivalence  ratio  as  a 
function  of  fuel  penetration  height  obtained  from  the  parametric  study.  Of 
particular  interest  is  the  result  for  near-zero  penetration  height,  t/r  ~  0. 
Under  this  condition,  for  a  global  equivalence  ratio  of  unity,  the  predicted 
recirculation  zone  fuel /air  ratio  is  extremely  fuel -rich,  with  4^  =  4.5. 

These  high  values  of  equivalence  ratio  produce  a  substantially  cooler  recircula¬ 
tion  region  than  would  be  obtained  for  conditions  nearer  stoichiometric,  and 
this  in  turn  leads  to  reduced  flame  propagation  rates  and  lower  overall 
combustion  efficiency  in  the  combustor.  Evaluation  of  the  effect  of  the  rich 
fuel/air  ratio  in  the  recirculation  zone  on  overall  performance  is,  however, 
complicated  by  difficulties  that  are  encountered  in  modeling  chemical  kinetics 
in  very  fuel -rich  regions. 

There  is  no  currently  available  finite-rate  chemical  kinetics  model  which 
can  provide  adequate  results  for  complex  hydrocarbon  fuels  under  fuel -rich 
conditions,  at  temperatures  at  which  the  effects  of  dissociation  manifest 
themselves.  The  modular  model  formulation  used  for  the  computations  described 
in  this  paper  includes  as  options  either  a  one-step  finite-rate  model  (Ref.  4) 
or  the  quasiglobal  kinetics  model  (Refs.  2,  14).  In  the  former  case,  the 
products  of  combustion  are  the  fully-reacted  species  Wfl  and  COg,  so  that  the 
effects  of  dissociation  are  ignored,  and  the  overall  reaction  rate  is  under¬ 
predicted.  The  quasiglobal  model  used  in  the  modular  formulation  is  valid 
for  complex  hydrocarbon  fuels  and  includes  the  effects  of  dissociation,  but 
is  not  appropriate  for  equivalence  ratios  much  above  stoichiometric.  Other 
models  have  been  proposed  for  fuel -rich  conditions,  such  as  the  four-step 
model  of  Dryer  and  Glassman,  Ref.  12,  but  this  model  ignores  dissociation 
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FIGURE  5.  Recirculation  Zone  Equivalence  Ratio  as  a  Function  of 
Penetration  Height 
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effects.  Thus  it  is  difficult  to  obtain  reliable  estimates  for  the  chemical 
kinetic  rates  that  pertain  in  combustors  in  which  significant  regions  of  the 
flowfield  are  highly  fuel-rich. 

This  difficulty  is  reflected  in  the  performance  results  shown  in  Fig.  6. 
While  the  trend  of  combustion  efficiency  as  a  function  of  fuel  penetration 
distance  obtained  through  the  use  of  the  one-step  global  kinetics  model  is 
as  expected,  the  magnitude  of  the  overall  change  is  relatively  small  and 
probably  incorrect.  The  reason  that  the  magnitude  of  these  results  can  be 
expected  to  be  incorrect  is  that  the  global  finite-rate  model  can  be  expected 
to  significantly  overpredict  the  recirculation  zone  temperature  for  fuel -rich 
conditions,  since  the  only  effect  of  fuel-rich  operation  that  this  model 
recognizes  is  a  dilution  of  the  products  of  combustion  with  excess  fuel  for 
<t>  »  1.  With  this  model,  the  reaction  rate  is  linear  with  respect  to  fuel 
concentration,  whereas  evidence  exists  that  for  fuel -rich  conditions  reaction 
rates  do  not  continue  to  increase  with  fuel  concentration,  and,  as  noted 
earlier,  dissociation  effects  are  also  ignored. 

Nevertheless  these  results  indicate  that  for  the  conditions  considered, 
low  values  of  fuel  penetration  distance  are  directly  coupled  to  reductions 
in  overall  combustor  performance.  Since  further  development  testing  using 
fuel  injection  schemes  designed  to  increase  fuel  penetration  resulted  in 
substantial  performance  improvements,  it  can  be  considered  that  the  trends 
demonstrated  by  this  modeling  effort  have  been  experimentally  verified. 

In  the  descriptions  of  the  modular  approach  to  the  development  of  combustor 
models  outlined  in  Refs.  1,  4,  6  and  7  it  has  been  stressed  that  one  of  the 
features  of  the  modular  approach  is  that  each  of  the  elements  involved  in  the 
overall  combustor  model  can  be  independently  developed.  Just  such  an  approach 
is  being  followed  at  SAI  where,  under  a  DOE-sponsored  program,  the  extension 
of  the  quasiglobal  model  to  fuel-rich  conditions  is  being  explored  using  the 
well -stirred  reactor  element  of  the  model.  One  result  of  this  work  is  shown 
in  Fig.  7,  on  which  predictions  using  several  different  kinetics  models 
are  compared  to  measured  temperature  data  obtained  in  a  laboratory  well- 
stirred  reactor  burning  toluene  under  both  fuel-lean  and  fuel-rich  conditions. 
All  of  these  results  are  for  the  same  stirred  reactor  residence  time. 

For  the  fuel-rich  condition.  Fig.  7  shows  that  the  extended  quasiglobal  model 
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FIGURE  6.  Predicted  overall  combustion  efficiency  as  a 
function  of  fuel  penetration  height. 


substantially  improves  the  stirred  reactor  temperature  prediction  as  compared 
to  either  the  original  quasiglobal  formulation  (Refs,  2,  11)  or  the  four-step 
model  (Ref.  12). 


Although  the  extended  quasiglobal  formulation  has  not  yet  been  fully 
developed,  an  indication  of  the  effect  of  using  this  approach  in  a  combustor 
calculation  is  demonstrated  by  Fig.  8,  which  shows  recirculation  zone  tempera¬ 
tures  for  the  conditions  of  the  preceding  section,  computed  using  the  four-step 
model  (Ref.  12),  the  one-step  model  (Ref.  4),  and  the  extended  quasiglobal 
formulation  (Ref.  13).  While  these  results  are  provisional  in  that  the 
extended  quasiglobal  model  has  not  at  this  point  been  verified  for  ramjet  fuels 
or  operating  conditions,  they  indicate  a  significant  effect  on  recirculation 
zone  temperature.  In  particular,  at  4>^z  =  3.0,  the  recirculation  zone  tempera¬ 
ture  predicted  using  the  extended  quasiglobal  model  is  some  300°K  lower  than 
that  obtained  from  the  one-step  model.  This  difference,  coupled  with  the 
radically  different  species  distributions  obtained  from  the  two  models,  can 
be  expected  to  exert  a  significant  influence  on  predicted  flame  propagation 
rates. 


In  summary,  a  comparison  of  modular  model  predictions  with  available 
data  shows  that  the  model  is  capable  of  predicting  both  the  entrainment  rate 
into  sudden-expansion  combustor  recirculation  regions  and  the  overall  fuel/air 
ratio  in  these  regions.  Both  of  these  parameters  are  crucial  to  investigations 
of  flame  stability  in  a  liquid  fueled  sudden-expansion  combustor.  An  appli¬ 
cation  of  the  approach  in  the  interpretation  of  ramjet  combustor  test  data 
has  also  been  demonstrated:  This  demonstration  shows  that  the  model  provides 
a  means  of  verifying  suspected  causes  of  combustor  performance  deficiencies. 

The  inadequacy  of  current  chemical  kinetics  models  under  very  fuel -rich  condi¬ 
tions  in  circumstances  under  which  the  effects  of  dissociation  are  apparent 
provides  a  current  limitation  on  the  use  of  the  modular  model  approach. 

However,  related  chemical  kinetic  model  development  work  which  can  remove  this 
limitation  has  been  described. 
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Combustor  recirculation  zone  temperatures 
computed  using  different  chemistry  models 


3.  MODIFICATIONS  TO  MODULAR  MODEL  FOR  BORON  SLURRY  FUEL  COMBUSTION 


Boron  slurry  fuels  have  a  high  energy  density  and  are  therefore  extremely 
attractive  for  volume- limited  ramjet  applications.  However,  the  advantages  of 
boron  slurry  fuels  are  only  available  if  combustion  efficiencies  comparable 
to  those  obtained  with  more  conventional  hydrocarbon  fuels  can  be  attained. 

In  past  combustor  development  work  the  realization  of  high  levels  of  combustion 
efficiency  has  been  achieved  only  at  the  cost  of  large  total  pressure  losses 
in  the  combustor,  which  forces  unacceptable  design  compromises. 

Because  the  mechanisms  involved  in  boron  slurry  combustion  are  considerably 
different  from  those  involved  in  the  combustion  of  liquid  hydrocarbon  fuels,  the 
empirical  design  guidelines  that  have  been  developed  for  conventional  fuels 
do  not  apply.  Even  for  conventional  fuels,  the  need  for  reduced  combustor 
volume  and  the  packaging  requirements  that  lead  to  the  use  of  a  sudden-expansion 
dump  combustor  have  together  reduced  the  applicability  of  empirical  design 
and  scaling  rules.  For  this  reason,  an  analytical  tool  for  the  design  and 
scaling  of  sudden-expansion  ramjet  combustors,  utilizing  conventional  hydro¬ 
carbon  fuels  has  been  developed  and  tested,  as  outlined  in  the  preceding 
section.  This  code  combines  a  semi -empirical  analysis  of  the  fuel  injection 
processes  with  a  sophisticated  turbulent  kinetic  energy  analysis  of  the  mixing 
process  and  a  quasiglobal  finite-rate  hydrocarbon  chemical  kinetics  mechanism 
in  an  iterative  analysis  of  a  sudden-expansion  combustor  with  large  scale 
embedded  recirculation  zones.  Because  of  current  interest  in  boron  slurry 
fuels  for  ramjets  and  the  inapplicability  of  empirical  guidelines  for  the 
design  of  ramjet  combustors  using  these  fuels,  the  development  of  a  further 
modification  to  the  modular  model  for  liquid-fueled  sudden-expansion  ramjets, 
to  extend  the  approach  to  boron  fuel  combustion,  has  been  initiated  under 
this  AFOSR  program. 

Boron  slurry  fuel  combustion  is  a  complex  phenomenon  involving  two-phase 
flow  and  finite-rate  reaction  processes  in  the  gas  phase  and  on  the  particle 
surface.  The  boron  is  generally  in  the  form  of  micron-size  particles  suspended 
in  a  carrier  hydrocarbon  fuel.  One  of  the  difficulties  involved  in  boron 
combustion  is  that  a  surface  oxide  layer  can  form  on  the  particle.  This  oxide 
has  a  boiling  point  below  that  of  the  metal,  so  that  the  boron  combustion 
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process  involves  a  diffusion  of  oxidizing  gases  through  the  oxide  layer  with 
a  surface  reaction  mechanism  involved  at  the  particle  surface.  Since  this  is 
a  relatively  slow  process,  one  strategy  for  obtaining  high  boron  combustion 
efficiencies  is  to  limit  the  amount  of  oxygen  reaching  the  surface  of  the 
particle  during  the  heat-up  of  the  particle,  i.e.,  carrying  out  the  particle 
heat-up  in  a  regime  which  is  fuel-rich  relative  to  the  hydrocarbon  carrier, 
and  then  exposing  the  hot  particles  and  remaining  fuel  to  additional  air. 

This  form  of  staged  combustion  is  difficult  to  achieve  in  practice,  and  the 
assembly  of  a  computer  analysis  to  aid  in  achieving  it  is  one  of  the  goals 
of  this  task. 

The  analysis  of  the  combustion  of  a  boron  slurry  fuel  in  a  ramjet  combustor 
begins  just  downstream  of  the  point  of  injection  of  the  fuel.  A  spray  is 
formed,  which  involves  a  distribution  of  droplet  sizes  within  which  there  is 
also  a  distribution  of  boron  particle  sizes.  As  the  fuel  evaporates  and  burns, 
there  remains  a  cloud  of  boron  particles,  with  varying  sizes  and  varying  oxide 
layer  thicknesses.  To  compute  the  subsequent  fate  of  these  particles  in  the 
combustor,  it  can  be  assumed  that  the  particles  are  in  near  dynamic  equilibrium, 
which  means  that,  with  respect  to  the  mean  velocity  in  the  combustor,  the 
particles  follow  the  gas  phase,  but  they  do  not  necessarily  follow  the  turbulent 
fluctuations  of  the  gas  phase.  Thus  there  is  an  effective  diffusion  of  parti¬ 
cles  relative  to  the  gas  phase. 

In  the  analysis,  the  continuous  distribution  of  particles  is  discretized 
into  a  number  of  classes,  with  the  class  described  by  an  average  particle 
size  and  oxide  layer  thickness  range.  It  is  further  assumed  that  these 
particles  form  a  dilute  continuum,  so  that  the  motion  of  each  class  can  be 
described  by  a  species  transport  equation  and  the  thermal  history  of  each 
class  by  a  particulate  energy  equation.  The  coupling  of  these  species  and 
energy  equations  to  the  gas  phase  combustion  and  heat  transfer  mechanisms  is 
through  source  terms  which  model  the  rate  of  consumption  and  heat  transfer 
phenomena  for  each  class  of  particles.  Since  the  boron  oxidation  process  is 
in  general  occurring  simultaneously  with  the  oxidation  of  the  hydrocarbon 
carrier  fuel,  these  chemical  kinetic  phenomena  must  be  considered  together, 
which  requires  the  addition  of  species  involved  in  the  boron  oxidation  process 
to  those  already  included  in  the  model. 
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The  two  major  modifications  to  the  modular  model  that  are  required  for 
the  modeling  of  boron  slurry  combustion  are  the  development  and  introduction 
of  a  particle-tracking  methodology  and  the  development  and  introduction  of 
submodels  for  the  boron  particle  heat-up  and  combustion  processes.  These  two 
modifications  are  the  subject  of  this  Section.  As  there  currently  exists  no 
comprehensive  model  for  boron  particle  combustion,  the  two  most  complete 
models  available  have  been  coded  as  interchangeable  subroutines  which,  in 
addition  to  their  use  within  the  modular  model,  can  also  be  utilized  as 
single-particle  models.  In  the  latter  mode,  these  subroutines  can  be  used 
in  the  development  of  more  complete  models  of  the  boron  combustion  process, 
through  comparison  of  single-particle  model  results  with  relevant  experimental 
data.'  The  particle  tracking  formalism  has  been  incorporated  into  the  code 
utilizing  a  total  of  ten  particle  categories:  as  each  additional  category 
involves  an  additional  species  transport  equation,  in  addition  to  an  increase 
in  the  overall  storage  required  for  the  code,  ten  categories  have  been  selected 
as  providing  the  best  balance  between  the  need  for  discrimination  of  particle 
size  effects  and  the  overall  code  run  time  required  for  solution. 

3.1  BORON  IGNITION  AND  COMBUSTION  MODELS 

There  are  a  number  of  models  for  the  ignition  and  combustion  of  single 
particles  of  boron.  A  variety  of  these  have  been  investigated  for  possible 
application  in  the  modular  model  for  slurry-fueled  ramjets,  including  those 
reported  by  King  (Refs.  14,  15),  Meese  (Ref.  16),  Vovchuck  (Ref.  17),  Mohan 
(Ref.  18)  and  Edelman  (Refs.  19,  20).  The  component  equations  used  in  each 
of  these  models  are  summarized  in  Table  1,  which  also  incorporates  the 
assumptions  built  into  the  models;  the  nomenclature  for  Table  1  is  given  in 
Table  2. 

In  reviewing  these  models,  the  emphasis  has  in  each  case  been  on  the 
utility  and  generality  of  the  approach  with  respect  to  the  question  of  addressing 
the  combustion  of  boron  particles  under  the  wide  range  of  conditions  appropriate 
to  ramjet  operation.  The  attention  given  to  both  diffusion  controlled  and 
kinetically  controlled  mechanisms  was  examined,  as  was  the  characterization 
(if  any)  of  oxide  removal  by  water. 
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In  general,  the  formulations  are  based  upon  the  conservation  equations  for 
mass  and  energy.  The  momentum  equation  is  not  explicitly  coupled  into  the 
equation  set  but  in  certain  of  the  models  the  relative  velocity  effect  is 
accounted  for  through  a  Nusselt  number  correction  on  heat  and  mass  transfer. 

In  most  of  the  formulations  the  energy  equation  is  the  central  element  of  the 
model,  with  the  size  history  and  other  transport  relationships  representing 
auxiliary  elements  of  the  model. 

King  considers  three  intervals  in  the  heat-up  of  the  particle:  below 
phase  change  of  the  boron,  at  phase  change  of  the  boron,  and  above  the  phase 
change  of  the  boron  (melting  point  2450°K).  Meese,  who  based  his  work  on  a 
1972  paper  of  King's,  also  includes  the  boil-off  of  the  oxide  layer.  For 
conduction.  King  assumes  the  particle  moves  at  the  velocity  of  the  surrounding 
gas,  whereas  Meese  allows-for  the  incorporation  in  his  calculation  of  the 
convective  heat-transfer  coefficient.  Vovchuk  includes  no  convective  heat 
transfer  and  in  Table  1,  this  section  for  Vovchuk  is  omitted.  This  is  the 
form  used  when  a  section  does  not  apply  to  an  author’s  work.  Edelman  also 
allows  for  convective  heat  transfer.  The  format  for  the  listing  of  the 
radiation  component  in  Table  1  is  the  same  as  that  used  for  the  convection  terms. 
Each  radiation  equation  could  be  interchanged  with  the  other  radiation  equations 
with  only  minor  modifications  to  match  an  author's  model  (i.e.  volume  based 
on  a  slab  or  a  particle).  It  should  be  noted  that  Vovchuk  includes  a  radiation 
term  which  he  eventually  neglects. 

For  oxide  formation,  Ki I®//  tfeese ,  andlToharTuse  the  global  reaction 
B  +  3/402~-*~  1/2B203.  Vovchuk  allows  for  lower  oxides  BO,  BO2,  B2O2,  forming 
at  the  surface.  Edelman  includes  BO  with  B203  but  if  BO  forms  it  does  not 
oxidize  further.  All  reactions  are  diffusion  controlled  for  the  first  four 
authors  and  kinetics  controls  the  reactions  Edelman  uses. 

For  oxide  consumption,  only  King  and  Meese  consider  the  chemical  effect 
of  water  on  boric  oxide  (8^0^  +  H20-*-2HB02) .  The  only  other  mechanism 
considered  for  oxide  removal  is  evaporation. 

Diffusion  equations  in  King's  work  are  algebraic  in  form  and  depend  upon 
absolute  values  of  pressures  at  infinity  and  the  surface.  Vovchuk  integrates 
the  partial  pressure  from  the  surface  to  infinity.  For  diffusion  through  the 
oxide  layer  King  and  Meese  have  similar  equations  that  vary  only  in  one  constant. 
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Meese's  constant  was  derived  from  data  unavailable  to  King  and  fits  the  data 
curve  better. 

A  major  part  of  King's  model  is  the  phase  change  of  the  boron.  Meese  has 
a  similar  equation  though  derived  in  a  slightly  different  fashion.  Edelman 
assumes  a  finite  evaporation  rate  of  the  boron.  For  evaporation  of  the  oxide 
King  and  Meese  are  again  similar.  Vovchuk  and  Edelman  have  no  oxide  layer 
to  remove.  Mohan's  phase  change  model,  instead  of  being  equilibrium  controlled, 
is  finite  rate. 

Edelman  is  the  only  author  who  includes  finite  rate  kinetics.  He  assumes 
the  particle  has  ignited  and  studies  surface  oxidation  kinetics.  Vovchuk  also 
studies  only  combustion  in  the  diffusion  controlled  mode. 

The  Nusselt  number  correlation  allows  for  the  effect  of  a  velocity 
difference  on  convective  heat  transfer.  King  sets  the  value  of  two  (Nu  =  2)  by 
assuming  a  zero  velocity  difference.  Meese  and  Edelman  allow  for  variation  in 
velocities. 

In  Meese's  model  nothing  else  occurs  while  the  oxide  layer  boils  off 
leading  to  an  independent  expression  for  the  oil  off  time.  The  end  product  of 
Vovchuk's  work  is  the  burning  time.  Mohan  also  derived  an  expression  for  the 
burning  time.  The  radius  and,  in  cases  where  the  oxide  is  considered,  the 
oxide  layer  thickness  equation  are  basic  to  all  but  Vovchuk's  model. 

Three  of  the  models  consider  ignition:  King,  Meese,  and  Mohan.  They  all 
essentially  take  it  to  be  the  point  where  the  oxide  layer  boils  or  flashes  off. 
Other  assumptions  and  restrictions  made  by  each  author  appear  at  the  end  of 
Table  1. 

Based  upon  this  review  and  formulation,  the  two  most  promising  models 
for  use  in  ramjet  applications  are  those  of  King  and  Edelman.  These  are  most 
general  in  their  respective  formulations;  the  others  are  either  derivatives 
of  them  or  are  too  restrictive  to  be  usefully  extended. 

3.2  BORON  PARTICLE  COMBUSTION  SUBROUTINES  FOR  THE  MODULAR  MODEL 

The  emphasis  in  this  phase  of  the  present  work  was  on  the  creation  of 
a  single-particle  combustion  code  which  could  employ  either  the  methodology 
of  King  (Refs.  14,  15)  or  that  employed  by  Edelman  (Refs.  19,  20).  The  concept 
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TABLE  1.  Summary  of  Available  Boron  Ignition  and 
Combustion  Models. 


(  J  *  rt  PS  fP8(S  or  t)  4  4  *  rP  1  PB£ Oj  CP«20j  )  "a/  *  RB  '  R£  <AHVAP^  '  *H  *  *  *  **>  *  X*’ 

•  [h  <r_  -  V  ♦  o  c  o,  (T^  -  T‘>] 

'"B  CPB  4  "B.Oj  CPB20jJ  W  *  *  °  '  RE  WVAP  4  4  'rJ  ‘  [h  ‘T9  -  V  4  0  »  ‘  <TW0  *  T?>] 


*02  4  ».  W  -  >.2  Vj(T)  ’  Q« 
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«*  - 1  ~ ' L 
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LOWER  OXIDES  ONLT  AT  SURFACE 
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Vovchuk 
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Mohan 

NO  OXIDE  LATER 
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PHASE  CHANGE 

BORON 
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D 
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PHASE  CHANGE 

OXIDE 
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D 
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NO  OXIDE  LATEX 

Vovcbuk 

K 

=-  ■  2  0  /d  C  a  LIMITED  DEPENDING  ON  a 

%  * 

Mo  Nan 

NO  OXIDE  LATE* 

Edc)m*n 

CHEMISTRY 

M2°(9)  4  B2°3U)  *  ?HB0Z(9>  »  °2(9>  *  * 

King 

H2°(9)  4  B2°3(l)  *  ZHB02(9)  *(*.»)  *  *  °2(«>  *  ?  BZ°3(i.9) 

Meet* 

?B02(9)~B2°3(9>  4  7°Z(9> 

B203(9)=  ZB0(9)  +i0Z(9) 

B2°2(9);;=ZB0(9) 

Vovchuk 
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Mohan 
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CHEMICAL  KINETICS 
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Me°(®>*  HB%>*  °z(g)*  b?°3U.9)*  Bt»,0*  "2(9) 


SPECIES  FOUND 
AT  THE  SURFACE 


SPECIES  FOUND 
OFF  SURFACE 
(GASES) 


H2°(g)*  HB°2(g)*  °2(g)*  B2°3(t.g)*  *(..»)•  "?(*) 


■(..I)*  “(g)*  *°2(g)‘  *2°2(g)*  B2°3{g)*  °2(g) 
B(0*  *2e2(g)'  °2(g)*  B2°3(i.g) 


B(*.0*  B2°2(g)’  °2(g)'  "2(g)*  B0(g) 


Edelman 


Kir*g 
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Ag  eq.  Includei:  an  Arrhenius  v tic os 1 ty> temperature  l«w,  k  *  A 

Henry's  la*.  Pt/XI  •  K1  , 

(♦  »„)S  T 

Wilke  correlation,  Ot,  •  7.4  •  10  *  - - - 

«  V  ••• 

nl  A 

King 

fQUAT ION  RtflRlHClS 

0  Stefan's  constant  •  5.669  •  10'*  W/«*  ®K  -  0.1714  •  10'*  BTU/h-ft*  °8* 

0  •  O0  (T/T0)‘  »*  A  -  Xe  (1/T#)*  « 

Vovchuk 

Rb  •  X  K  0  /  h  Pe  -  exp  (H/A*T#  -  H/*#TC> 

K  0  can  be  expressed  1ft  an  Arrhenius  form  (see  King) 

d  In  p  l 

Clausius  -  Claperyon  — -- — -  • 

at  HI 

BBS 

IQUATIONS  FOR  h 

h  -  0.347  •  10'*  (Nu)  I*  */  (rf  *  X)  (1974)  h  -  0.694  •  10'‘  !•■•/  (1972) 

King 

h  -  (1  ♦  0.296  Pr1/’  Re’/*)  k  /r 

9 

Heese 

The  first  stage  of  a  t*o  stage  burn.  During  this  time  the  oxide  layer  Is  removed. 

King 

IGNITION 

A  state  at  which  no  energy  barriers  stand  between  that  state  and  steady-state  combustion. 

I 

Only  considers  combustion. 

mm  I 

Similar  to  King's. 

Mohan 

N/A 

fcfel  mu 

Oxide  generation  controlled  by  diffusion  of  oxygen  through  liquid  oxide  layer. 

Uniform  particle  temperature. 

Vaporisation  and  diffusion  considered  as  series  resistances. 

H^O  removes  B2O3  by  diffusion  limited  reaction  (endothermic). 

Heat  gain/loss  by  convection,  radiation. 

Partial  pressure  of  HB02  negligible  in  the  free  stream. 

Particle  is  stationary  with  respect  to  surrounding  gasses. 

Paper  only  deals  with  ignition. 

King 

ASSUMPTIONS 

Initial  B2O3  layer  thickness  is  uniform. 

Uniform  particle  temperature. 

The  convective  heat  transfer  coefficient,  h,  is  known. 

Surroundings  treated  as  black  at  temperature  Traq. 

Prior  to  oxide  removal,  boron  reaction  rate  is  diffusion  controlled  and  forms  62O3. 

B2O3  evaporation  follows  an  Arrhenius  type  relation. 

Particle  temperature  remains  constant  during  melting  and  boiling. 

Variation  in  h,  surface  reaction  rate  and  oxide  evaporation  rate  are  neglected  during  oxide 
boiling. 

Water  vapor  reduces  the  activation  energy  of  oxide  evaporation. 

Initial  particle  temperature  is  below  oxide  boiling  and  boron  melting. 

B2O3  formed  during  boiling  may  be  neglected. 

No  net  radiative  exchange  between  particle  and  cloud. 

Convective  heat  transfer  neglected  In  final  computations. 

Wees* 

The  particle  has  been  heated  to  equilibrium  temperature. 

Reaction  occurs  at  the  surface. 

The  diffusion  coefficients  for  the  oxidant  and  reaction  products  are  identical. 

Radiation  can  be  neglected  to  a  first  approximation  (d  «  100pm). 

Particle  temperature  Is  above  oxide  boiling  and  boron  melting  points. 

Vovchuk 

Diffusion  controlled  B  ♦  3/4  Oj  ■»  177  B^Oj. 
liquid  layer  is  Isothermal. 

Ideal  gas. 

Absorptive  equilibrium  for  Oj  at  C  and  ID  steady-stage  diffusion  of  0?  within  liquid  layer. 

Plane  two-dimensional  slab  model. 

Mohan 

No  oxide  layer. 

Uniform  particle  temperature. 

Surface  combustion  kinetics. 

Mnite>rate  evaporation  of  boron  from  the  surface. 

Cdelsan 
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TABLE  2.  Nomenclature  for  Expressions  Summarized  in  Table  1 


2 . 1  King's  Model 


'PB1 

'PBS 

'pb2o2 


liquid  boron  heat  capacity  (cal/gm°K) 
solid  boron  heat  capacity  (cal/gm°K) 
liquid  boron  oxide  heat  capacity  (cal/gm°K) 


f 

h 

k 


(mw)8 

(mw)b2o2 

Nu 


P 

Pr 


b20j,  surf 
P° 

B2°3 

PH20,  surf 
PHB02,  surf 
- 

qrx 

QRX2 


fraction  of  boron  in  the  liquid  phase 

gas-particle  heat  transfer  coefficient  (cal/cm  sec°K) 

mass  transfer  coefficient  for  transport  of  boric  oxide  gas  from 

p 

particle  to  free-stream  (gm-mol/cm  atm  sec) 

boron  atomic  weight  (gm/gm-mol) 

boric  oxide  molecular  weight  (gm/gm-mol) 

Nusselt  Number 
total  pressure  (atm) 

oxygen  partial  pressure  in  free  stream  (atm) 

boric  oxide  partial  pressure  adjacent  to  particle  surface  (atm) 

boric  oxide  vapor  pressure  (atm) 

water  gas  partial  pressure  adjacent  to  particle  surface  (atm) 
HB02partial  pressure  adjacent  to  particle  surface  (atm) 
boric  oxide  partial  pressure  in  free-stream  (atm) 
heat  release  of  B(s)  +  |o2  -*•  ^O^O)  (cal/gm-mol) 
heat  release  of  6(1)  +  |o2  -*■  ^BgO^O)  (cal/gm-mol) 

3 

gas  law  constant,  82.06  atm-cm  /gm-mol°K 

molar  rate  of  boron  consumption  (gm-mol/sec) 

molar  evaporation  rate  of  boric  oxide  (gm-mol/sec) 

molar  rate  of  removal  of  by  water  reaction  (gm-mol/sec) 

boron  particle  radius  (cm) 

free  stream  gas  temperature  (°K) 

particle  center  temperature  (°K) 

particle  temperature  (°K) 

surroundings  radiation  temperature  (°K) 
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I 


^2^3  *^2 


uhbo2,n2 

dh2o,n2 

AH., 


pb2o3 


particle  surface  temperature  (0K) 
molecular  volume  of  species  j  (cm  ) 
oxide  layer  thickness  (cm) 

o 

diffusivity  of  gaseous  boric  oxide  in  nitrogen  (cm  /sec) 

2 

diffusivity  of  HB02  in  nitrogen  (cm  /sec) 

2 

diffusivity  of  water  gas  in  nitrogen  (cm  /sec) 

heat  absorbed  by  reaction  of  H20  with  BgO^O ) (cal/gm-mol ) 
heat  of  fusion  of  boron  (cal/gm) 
heat  of  vaporization  of  B^^OMcal/gm-mol) 
evaporation  coefficient  of  boric  oxide  liquid 
surroundings  absorptivity 

2 

thermal  diffusivity  of  boron  (cm  /sec) 

1?  2  4 

Stefan-Bol tzmann  constant,  1.354*10  cal/cm  sec°K 
particle  emissivity 

Hertz-Knudsen  impingement  factor  (gm-mol/cm  atm  sec) 
time  (sec) 

boron  density  (gm/crn  ) 

3 

boric  oxide  density  (gm/cm  ) 


2 . 2  Model  of  Meese  and  Skifstad 

2 

A  surface  area,  cm 


specific  heat,  cal/g-°K 

particle  diameter,  cm 

activation  energy  for  evaporation 

heat  of  fusion  of  boron,  cal/g 

heat  of  vaporization  of  boron  oxide,  cal /mole 

convective  heat-transfer  coefficient,  cal/cm  -sec-°K 

thermal  conductivity,  cal/cm-sec-°K 

molecular  weight,  g/g-mole 

mass,  g 

pressure,  atm 

Prandtl  number  of  the  gas 

heat  of  combustion  of  boron  with  diatomic  oxygen  to  form  liquid 


BgOj,  cal/mole;  total  energy,  cal 
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R 

re 

Re 

r 

Ar 

T 

t 

xh2o 

a 

e 

X 

P 

a 


molar  rate  of  consumption  of  boron,  g-mole/sec;  gas  constant 
molar  rate  of  evaporation  of  boron  oxide,  gm-mole/sec 
Reynolds  number  based  on  particle  diameter  and  relative  velocity 
between  the  gas  and  particle 
radius,  cm 

oxide  film  thickness,  r2-r.| ,  cm 
temperature,  °K 
time,  sec 

mole  fraction  of  present  in  the  gas 

surrounding  absorptivity  (=1.0) 
particle  emissivity  (=0.8) 
particle  burning  rate,  cm/sec 
density,  g/cm 

2  4 

Stefan-Boltzmann  constant,  cal/cm  sec°K 


2.3  Mohan's  Model 


C 

Ct 

D 

d 

h 

k 

k' 

L 


average  molecular  velocity  normal  to  the  surface  (cm/sec) 

heat  capacity  of  B203(L)(cal/gm°K) 

diffusion  coefficient  (cm^/sec) 

diameter  of  particle  (cm) 

thickness  (cm) 

effective  distribution  coefficient  related  to  the  ratio  of  the 

mole  fraction  of  02  absorbed  in  B203(Jt)  to  gas  phase  mole  fraction 

distribution  ratio  (liquid  phase  to  gas  phase  0,  concentration) 

sum  of  heat  losses  by  radiation  and  conduction  — 2 - 

cal  sec 


P 

qB 

qC 

^B  or  C 
R° 

T 

t 

W1 

X 


total  ambient  pressure  (atm) 

heat  liberated  per  unit  mass  of  boron  consumed  at  B  (cal/g) 

heat  absorbed  per  unit  mass  of  B203  vaporized  at  C  (cal/g) 

regression  rate  of  the  surface  (cm/sec) 

universal  gas  constant  (8.316*10"^erg/gm-mole  °K) 

temperature  (°K) 

time  (sec) 

atomic  weight 

mole  fraction  of  02  in  the  ambient  atmosphere 

ratio  equilibrium  vapor  pressure  of  B203  at  Tc  to  total  pressure 
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a  evaporation  coefficient 

e  emissivity 

A  gas  phase  thermal  conductivity  (cal/cnrsec  °K) 

Pg  q  density  boron  and  (gm/cm) 

x  time  for  the  new  combustion  regime  to  appear  (sec) 

B,C  surfaces 


2.4 


do 


AH. 


M„ 


Pi 


R 

r 

T 

€ 

A 


Vi 


a 


T 


r 


Vovchuk's  Model 

molar  specific  heat 

diffusion  coefficient 

initial  diameter  of  the  boron  droplet 

total  enthalpy  of  species  i 

heat  of  reaction  of  species  i 

molar  flux  of  species  i 

mass  flow  of  the  boron  reacting  in  unit  time 
total  pressure 

partial  pressure  of  species  i 

radiative  heat  flux 

universal  gas  constant 

radius 

temperature 

emissivity 

thermal  conductivity  of  gas 
molecular  weight 
stoichiometric  coefficients 
density 

Stefan's  constant 
droplet  combustion  time 


Subscripts 

B  (molten)  boron 

K  particle  or  particle  surface 

•  infinity 


Note:  units  not  stated,  therefore  any  consistent  set  of  units  may  be  used. 
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number  of  moles  of  oxygen  per  mole  of  product 
the  total  (dimensionless)  particle  consumption  rate 
number  of  moles  of  boron  per  mole  of  product 
drag  coefficient 

specific  heat  of  the  gas  (ft^/sec^  °K) 

2  2 

specific  heat  of  boron  at  the  surface  (ft  /sec  °K) 

heat  of  reaction  (calorie/gram) 
surface  reaction  rate  constant 
latent  heat  of  vaporization 
molecular  weight  of  boron  (g/g-mole) 
molecular  weight  of  nitrogen  (g/g*mole) 
molecular  weight  of  oxygen  (g/g-mole) 

molecular  weight  of  the  product  (g/g-mole) 
average  molecular  weight  at  the  particle  surface  (g/g-mole) 
pressure  (in  atmospheres) 
vapor  pressure  of  boron  (atm) 

particle  radius  (feet) 
gas  temperature  (°K) 
particle  temperature  (°K) 
temperature  of  the  wall  (°K) 
gas  velocity  (ft/sec) 
particle  velocity  (ft/sec) 

mass  fraction  of  boron  at  the  surface  of  the  particle 
mass  fraction  of  nitrogen  at  the  surface 

mass  fraction  of  oxygen  at  the  surface  of  the  particle 

mass  fraction  of  oxygen  in  the  surrounding  (gas) 

mass  fraction  of  the  product  at  the  particle  surface 
fraction  of  vaporized  product 
particle  bulk  density  (slugs/ft  ) 
emissi vity 

the  thermal  conductivity  of  nitrogen  (lbs  force/sec  °K) 

o 

gas  density  (slugs/ft  ) 
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p  density  at  the  surface  of  the  particle  (slugs/ft  ) 

j  » 

o  Stefan  Boltzmann  constant  (lbs  force/ft*sec  °IC) 
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was  to  develop  a  code  which  could  stand  alone  as  a  tool  to  be  used  in  the 
investigation  of  single-particle  combustion  phenomena,  but  which  could  also 
be  easily  incorporated  into  the  modular  model  format  without  substantial 
modification.  This  dual  capability  was  achieved  by  structuring  the  code  as 
a  series  of  subroutines  which  do  not  depend  on  the  basic  code  in  which  they  are 
found.  In  the  stand-alone  code  they  are  called  out  of  a  main  routine  which 
does  little  more  than  obtain  the  needed  inputs  and  call  the  first  of  them;  in 
the  modular  approach  this  calling  subprogram  is  a  simple  routine  which  corrects 
existing  data  to  the  units  and  variables  required  and  calls  the  boron  particle 
subroutine  package.  The  subroutine  package  itself  needs  no  modification  in 
either  case. 


Both  models  involve  the  solution  of  a  set  of  coupled  ordinary  differential 
equations.  The  integration  was  carried  out  using  subroutine  HPCG,  one  of  the 
standard  differential  equation  solvers.  HPCG  uses  Hamming's  modified  predictor- 
corrector  method  to  obtain  an  approximate  solution  to  the  system.  Although 
not  self-starting  (the  initial  solution  is  obtained  using  a  Runge-Kutta 
technique)  the  technique  is  superior  to  the  Runge-Kutta  in  that  it  requires 
evaluation  of  the  right-hand  side  only  twice  per  step.  Furthermore,  because 
the  local  truncation  error  falls  out  automatically  calculation  of  the  proper 
step  size  is  facilitated  (it  has  the  capability  of  moving  in  either  direction 
as  the  solution  proceeds,  but  it  never  exceeds  an  input  maximum).  The 
technique  used  is  a  stable  fourth  order  scheme. 

3.2.1  King's  Model 

3. 2. 1.1  Model  Description 

King's  model  includes  the  inhibiting  effect  of  the  boric  oxide  layer  on 
ignition  and  explicitly  calculates  the  particle  melt  time  before  ignition. 

It  is  thus  the  better  model  of  the  ignition  process.  However,  because  it 
does  not  include  the  kinetically  controlled  terms,  it  is  poor  for  predicting 
results  for  ignited  particles  (indeed,  it  predicts  a  total  thermal  runaway). 

The  basic  model  consists  of  four  coupled  equations  relating  particle 
temperature,  particle  size,  oxide  layer  thickness  and  the  fraction  of  the 
particle  which  is  melted.  It  can  be  expressed  as  follows: 


*BQRX  ~  REAHvap  "  RHaHH  +  47,(rp+x)2  ^h(VTp)+qcaR^Trad~Tp^  T  <245Q  f=( 


4  n  r„  Pn  C  +  4ur2X  pR  n  C 
3  p  PB  pBS  p  B,0,  p 
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^  = 
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0,  Tp  =  2450,  f<l 
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QRXZ 

rb 

rh 
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rad 
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aR 

AHH 


liquid  boric  oxide  heat  capacity  (cal/gm  °K) 

liquid  boron  heat  capacity  (cal/gm  °K) 

solid  boron  heat  capacity  (cal/gm  °K) 

gas-particle  heat  transfer  coefficient  (cal/gm  °K) 

3  1 

heat  release  of  B(S)  +  J  ^  2  B2®3^^ca^9m_m0^e^ 

heat  release  of  B(L)  +  f  ^  1  ^ (c:al/gm-mo7e) 

molar  rate  of  boron  consumption  (gm-mol e/sec) 

molar  rate  of  removal  of  B203  by  water  reaction  (gm-mole/sec) 

particle  radius  (cm) 

particle  temperature  (°K) 

surroundings  radiation  temperature  (°K) 

free  stream  gas  temperature  (°K) 

oxide  layer  thickness  (cm) 
surroundings  absorptivity 

heat  absorbed  by  the  reaction  of  H20  with  B203U) (cal /gm-mol e) 
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AHvap  heat  of  vaporization  of  (cal/gm-mole) 

c  particle  emissivity 

O 

pg  boron  density  (gm/cm  ) 

3 

pR  n  boric  oxide  density  (gm/cm  ) 

b2u3 

o  Stefan-Boltzman  Constant  (cal/cm^sec  °K4) 

0,  Tp<2450,  f=0 

RB^RX~^EAHvap~RHA^H+4lT ^  Vx ^ 2^h ^ T°°~Tp^+aEaR^Trad~Tp \ 2 ,  Tp  =  2450,  f<l 

4/3^pPBAHm 

VO,  Tp>2450,  f=l  (4) 

where  AHm  =  heat  of  fusion  of  boron  (cal/gm),  and  the  other  variables  are  as 
defined  above. 


^2.  = 
dt 


Rb(mw) 


B 


4"rpcB 


(5) 


where  (MW)g  =  boron  atomic  weight  (gm/gm-mole) ,  and  the  other  parameters  are 
as  above.  Finally, 


dx 

dt 


(Rb/2  -  Re  -  Rh)(MW)B2Q3 

4,rpCB203 


(6) 


where  (MW)D  n  is  the  boric  oxide  molecular  weight  (gm/gm-mole)  and  all  other 
h2u3 

variables  are  already  defined. 

As  can  be  seen,  the  King  model  is  actually  three  models,  of  which  one 
is  selected  based  on  current  particle  temperature,  and  melt  history.  This 
requires  careful  handling  in  the  coding.  The  HPCG  routine  needs  a  subroutine 
to  evaluate  the  derivatives  at  each  time  step.  One  cannot,  however,  merely 
wait  until  the  particle  temperature  meets  some  criterion,  then  start  using  the 
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melt  equation  because  the  integrator  may  call  this  routine  many  times  for  every 
final  solution  it  accepts.  Furthermore,  there  is  no  easy  way  to  tell  while  in 
the  derivative  routine  whether  or  not  the  solution  this  time  will  be  accepted. 
Thus  it  was  necessary  to  write  three  essentially  separate  routines  for  the 
three  states  of  the  particle  heating,  melt,  and  burn. 

Initially  the  time  range  of  interest  is  passed  to  the  integrator,  which 
starts  with  the  assumption  that  the  particle  is  not  melted.  When  (and  if) 
the  particle  reaches  (or  passes)  the  melt  temperature  (2450°K)  the  integration 
is  stopped,  values  are  set  back  to  those  at  which  Tp  =  2450  exactly,  and  the 
integration  procedure  is  initialized  for  the  new  problem  (defined  now  as  the 
calculation  of  the  fraction  of  the  boron  which  is  liquid,  f  (Eq.  4)).  When 
the  liquid  fraction  reaches  (or  passes)  unity,  the  integration  is  again  stopped 
and  the  values  are  again  reset  (this  time  to  those  at  which  the  liquid  fraction 
is  exactly  1.0)  and  the  integration  is  initialized  for  the  burn  problem. 

This  methodology  allows  automatic  correct  handling  of  those  cases  in  which 
the  particle  never  melts,  those  in  which  only  partial  melting  occurs,  as  well 
as  those  cases  in  which  ignition  occurs.  The  implementation  is  all  automatic 
and  transparent  to  the  user. 

Because  the  individual  time  steps  are  small  with  respect  to  the  variables 
being  calculated  it  was  decided  that  a  linear  interpolation  gave  sufficient 
accuracy,  and  represented  the  best  approach  to  the  problem  of  resetting  the 
variables  at  the  boundaries  of  the  particle  states.  A  comparison  of  King's 
results  with  those  obtained  using  the  SAI  code  verifies  the  adequacy  of  this 
approach. 

The  approach  was  checked  by  recalculating  some  of  King's  original  results 
using  the  King  model  code.  Two  cases  were  selected  for  this  on  the  basis  of 
near-critical  temperature  determining  whether  or  not  ignition  occurred. 

In  both  causes  the  particle  radius  is  10  microns  and  the  oxide  layer  is  0.1 
microns.  The  initial  particle  temperature  is  1800°K;  the  oxygen  mole  fraction 
is  0.2  and  the  pressure  is  5  atmospheres. 

The  first  case  has  a  gas  temperature  of  2100°K.  Fig.  9  shows  the 
particle  temperature  histories  as  calculated  by  King  and  by  SAI.  The  results 
are  essentially  identical.  The  observed  differences  are  more  likely  due  to 
inaccuracies  in  reading  the  King  curve  from  Ref.  16  than  to  be  actual  numerical 
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FIGURE  9.  Temperature  history  for  the  2100°  case. 


differences.  Fig.  10  compares  the  two  calculations  of  the  oxide  layer  thickness 
as  a  function  of  time.  Again  it  is  seen  that  the  two  calculations  are  essentially 
identical.  Fig.  11  gives  the  particle  radius  history  for  the  SAI  calculation 
(King  did  not  present  this  result). 

The  second  comparison  case  was  run  for  an  identical  case  except  that  the 
gas  temperature  was  2000°K.  Here  King  shows  that  true  ignition  does  not  occur. 
Fig.  12  shows  that  the  SAI  code  predicts  essentially  the  same  result  as  was 
obtained  by  King.  Fig.  13  also  confirms  this  result  for  the  oxide  layer 
thickness  (the  apparent  difference  again  may  well  be  due  to  difficulty  in 
reading  the  detail  on  the  plot  of  King's  results  -  if  the  data  were  plotted 
in  the  same  scale  as  they  appear  in  Ref.  16  there  would  be  no  discernable 
difference).  Comparison  of  Figs.  11  and  14  shows  very  little  difference 
initially.  However,  in  the  degenerate  combustion  case  the  initial  trend  is 
continued  whereas  in  the  full  ignition  case  (Fig.  11)  the  effect  of  particle 
heatup  and  melting  is  evident. 

The  two  cases  taken  together  demonstrate  that  the  King  model  has  been 
faithfully  reproduced  in  the  SAI  code.  The  system  behaves  just  as  predicted 
by  King  on  both  sides  of  the  critical  gas  temperatures. 

3.2.2  Edelman's  Model 

Edelman's  model  (Refs.  20,  21)  neglects  the  inhibiting  effect  of  the 
boric  oxide  layer  on  ignition  and  fails  to  calculate  the  particle  melt  time. 
However,  it  includes  all  the  kinetically  controlled  mechanisms  in  the  combustion 
equations.  Consequently,  it  is  superior  to  King's  model  for  examination  of 
what  happens  after  the  particle  has  ignited. 

The  full  Edelman  model  consists  of  nine  coupled  ordinary  differential 
equations.  However,  as  expanded  by  Edelman  (Ref.  21)  the  system  can  be 
simplified  to  three  equations  without  significant  loss  of  accuracy;  the 
remaining  equations,  although  actually  functions  of  time,  can  be  treated 
algebraically.  It  is  the  simplified  version  which  is  here  presented.  The 
three  differential  equations  are: 


ft  .  "wvv* 

dt  '  8R6$ 
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FIGURE  10.  Oxide  layer  history  for  the  2100°  case. 
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FIGURE  12.  Particle  temperature  history  in  the  2000°  case. 
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FIGURE  13.  Oxide  layer  thickness  history  for  the  2000°  case. 


FIGURE  14.  Particle  radius  history  in  the  2000°  case. 
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where 


and 


dR 

dt 


where 


5>  . 

dt 


where 


Cq  =  drag  coefficient 

Pg  =  gas  density  (slugs  per  cubic  foot) 

Vg  =  gas  velocity  (ft/sec) 

Vp  =  particle  velocity  (ft/sec) 

R  =  particle  radius  (feet) 

=  particle  bulk  density  (slugs/ft  ) 

-axH 

2R6~C~’  ^ 

s  p 

a  =  the  total  (dimensionless)  particle  consumption  rate  (see  Eq.  10) 
x  =  the  thermal  conductivity  of  nitrogen  (in  pounds  force/sec  °K) 

Nu  =  Nusselt  number 

2  2 

Cp  =  specific  heat  of  the  gas  (in  ft  /sec  °K),  and  the  remaining 
variables  as  already  defined. 
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specific  heat  of  boron  at  the  surface  (ft  /sec  °K) 

gas  temperature  (°K) 

particle  temperature  (°K) 

latent  heat  of  vaporization 

number  of  moles  of  oxygen  per  mole  of  product 

surface  reaction  rate  constant 

o 

density  at  the  surface  of  the  particle  (slugs/ft  ) 
mass  fraction  of  oxygen  at  the  surface  of  the  particle 

molecular  weight  of  boron  (g/g*mole) 

heat  of  reaction  (calorie/gram) 
molecular  weight  of  oxygen  (g/g-mole) 
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o  =  Stefan-Boltzmann  Constant  (lb* force/ft* sec  °KH) 
e  =  emissivity 

T  =  temperature  of  the  wall  (°K) 

VV 

and  other  parameters  are  as  defined  above. 

The  values  for  the  remainder  of  the  time-varying  parameters  are  treated  as 
simple  algebraic  equations  (using  the  current  value  of  Tp  and  R)  as  follows 


88.662C  aR  /mT  M, 

_ 1 _  p  _  JL  v  p 

ANu  \T  V*B(T  )  Mb  YB,S  k 


2MBKsV  cCpPsRA 


Mn  ANu 
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where  a  =  fraction  of  vaporized  product 

Py  j  =  vapor  pressure  of  boron  (atm) 

=  average  molecular  weight  at  the  particle  surface  (g/g-mole) 
Yg  s  =  mass  fraction  of  boron  at  the  surface  of  the  particle 
P  =  pressure  (in  atmospheres) 

and  other  variables  are  as  previously  defined. 
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where  Yn  =  mass  fraction  of  oxygen  in  the  surrounding  (gas) 
»°° 

B  =  number  of  moles  of  boron  per  mole  of  product 
and  the  remaining  variables  are  as  previously  defined. 
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where 


mass  fraction  of  the  product  at  the  particle  surface 
molecular  weight  of  the  product  g/g*mole) 


Y 

prod.s 
^prod 

and  the  remaining  parameters  are  defined  above 


=  average  molecular  weight  at  the  surface  = 


+ 


Y 

prod,s  J 
MProd 


where  =  molecular  weight  of  nitrogen  (g/g-mole) 

n2 

Yn  =  mass  fraction  of  nitrogen  at  the  surface 


and  the  remaining  variables  as  previously  defined. 
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The  Edelman  model  is,  in  a  sense,  more  easily  implemented  than  King's 
because  it  does  not  treat  separate  realms  of  phenomena.  However,  as  will  be 
seen  from  an  examination  of  Eqs.  10-13,  virtually  everything  is  a  function  of 
a,  which  itself  is  a  function  of  both  Yn  _  and  Y0  ...  Hence  one  is  faced  with  a 

Uo  »S  D  >  5 

set  of  three  highly  non-linear  equations^which  are  strongly  coupled. 

In  the  original  work  done  by  Edelman  an  iterative  solution  technique  was 
implemented  which  depended  on  the  ratio  of  surface  mass  fractions  to  surround¬ 
ings  mass  fractions  to  converge  to  a  solution.  However,  all  attempts  to 
implement  this  technique  were  met  with  uniform  failure.  Finally  this  approach 
was  abandoned  and  a  Newton-Raphson  scheme  was  derived  for  the  solution  of 
the  set  of  equations.  This  approach  quickly  led  to  convergence.  However,  a 
few  troublesome  anomalies  appeared,  the  most  serious  of  which  was  a  calculated 

negative  mass  fraction  for  the  product  in  some  areas.  Since  this  was  extremely 
-9 

small  (of  order  10  )  it  was  decided  to  be  due  to  inherent  inaccuracy  in  the 

numerics  (and/or  machine  roundoff)  and  zeroed  out  when  it  occurred  in  a  final 
converged  result  (but  not  before).  Consequently,  the  final  convergence 
criteria  were:  a  change  in  all  three  variables  from  one  step  to  another  of 
less  than  10'^  or  a  change  in  all  three  variables  of  less  than  0.1%.  If  either 


of  these  criteria  were  met,  the  solution  was  considered  to  have  converged 

unless  a  negative  mass  fraction  occurred.  A  significant  negative  mass  fraction 

causes  the  routine  to  iterate  again.  (A  significant  negative  is  defined  to 
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be  greater  in  absolute  value  than  10  ). 

Initial  results  from  the  Edelman  model  suggest  that  the  diffusion  terms 
are  not  the  dominant  ones  in  the  equations.  Hence  the  results  suggest  that 
the  driver  mechanism  is  kinetics;  indeed.  Fig.  15  shows  that  the  particle 
radius  decline  is  very  close  to  linear  with  time.  Even  the  temperature  curve 
for  the  J7.3  micron  particle  is  reasonably  close  to  linear  (Fig.  16).  Results 
for  a  larger  particle  (50  micron  diameter)  are  also  nearly  linear  (Fig.  17), 
and  a  plot  of  the  square  of  particle  diameter  vs.  time  (Fig.  18)  does  show  some 
curvature. 

3.3  PARTICLE  TRACKING  FOR  SLURRY  FUEL  COMBUSTION  ANALYSES 

The  combined  effects  of  continuous  reaction  of  boron  particles  and  the 
non-uniform  velocity  and  temperature  fields  within  the  combustor  provide 
conditions  resulting  in  a  spectrum  of  particle  sizes  at  every  point  in  the 
flow.  To  completely  specify  the  particle  reaction  rate  and  the  production  of 
gas  phase  species,  the  local  particle  size  and  number  density  distribution 
must  be  known.  This  in  turn  requires  the  development  of  a  methodology  for 
tracking  the  particle  history  throughout  the  flow. 

To  track  the  particle  history  through  the  flow,  the  total  amount  of  parti¬ 
culate  material  present  at  every  point  in  the  field  is  divided  into  a  set  of 
categories.  These  categories  can  be  based  on  any  parameter  which  can  be  used 
to  differentiate  particles:  for  the  boron  application  both  particle  size 
and  the  thickness  of  the  liquid  layer  (X)  on  the  particle  surface  are 
relevant  parameters.  Both  parameters  may  be  used,  with  liquid  layer  thickness 
forming  a  subclassification  in  a  basic  particle  size  classification,  but  it 
should  be  kept  in  mind  that  each  category  used  introduces  a  new  variable  in 
the  overall  combustor  model  with  associated  transport  equation  and  storage 
requirements. 

For  simplicity,  consider  a  categorization  based  on  particle  size.  Within 
each  size  range,  all  particles  are  combined  to  form  an  average  with  a  volume 
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FIGURE  IS.  Particle  radius  history  for  the  17.3  p  particle 
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INITIAL  CONDITIONS:  T„  =  2460  K,  P  !  1  atm 


FIGURE  17.  Particle  radius  history  for  the  25  y  particle. 


of  particle  diameter  as  a  function  of  time. 


to  surface  area  ratio  equal  to  that  of  the  summation  over  all  contributions. 
This  equivalent  particle  size  is  the  Sauter  mean  value  and  is  expressed  as 


.  £n 


,rf3/£V(2 


(16) 


where  n^  particles  have  radius  r . .  The  mass  fraction  of  each  category  is 
equal  to  the  sum  of  the  mass  fractions  of  each  contribution.  With  the  average 
size  of  category  particles  and  the  total  category  mass  fraction  known,  the 
appropriate  number  of  equivalent  particles  is  established  for  every  class. 

This  representation  provides  an  exact  equivalence  for  the  subsequent  rate 
process  provided  the  rate  of  change  of  mass  depends  directly  upon  the  particle 
surface  area. 

The  particle  distribution  array  in  terms  of  size,  mass  fraction,  and 
number  in  each  category  at  every  grid  point  can  only  be  altered  by  changes 
in  size  of  the  particle  due  to  chemical  reaction  (which  can  also,  for  boron, 
involve  changes  in  liquid  layer  thickness)  and  by  diffusion  of  particles  from 
the  neighboring  flow.  The  radius  of  the  particle  changes  in  accordance  with 
the  particle  consumption  expression  as  the  flowfield  develops.  Thus,  after 
every  step,  the  new  particle  size  is  given  by 


r  =  r  + 
o 


(17) 


Where  the  step  size.  Ax,  is  evaluated  based  on  the  stability  criteria  incor¬ 
porated  in  the  modular  model.  After  every  step  the  updated  radius  array 
must  be  reclassified  in  order  to  account  for  the  transition  of  particles  from 
one  category  to  another.  Since  it  is  possible  for  more  than  one  contribution 
to  be  made  to  a  particular  category,  the  total  category  mass  must  then  be 
reaveraged  using  the  Sauter  criteria,  to  establish  the  proper  radius  and  number 
for  the  group.  The  particle  array  in  terms  of  mass  fraction,  size  and  number 
density  resulting  from  the  kinetics  of  the  combustion  process  will  be  defined  as 

G(J,I)k 

R(J,I)k  at  grid  point 
N(J,I)k 


Each  category  in  the  particle  array  is  altered  by  the  diffusion  of  particles 
of  the  same  category  from  neighboring  upstream  locations.  Subdividing  the 
particle  mass  fraction  a.  into  each  separate  category  and  treating  the  category 
mass  fractions  as  separate  species,  the  diffusion  equation  can  be  written  in 
finite  difference  form  in  von  Mises  coordinates  as 


ai ,n+l  ,m 


,n,m 


+ 


Ax 

/(A*)2 


(18) 


where  the  n,m  subscripts  refer  to  the  grid  point  locations  shown  in  Fig.  19. 


The  second  term  in  Eq.  18  represents  the  change  in  due  to  the  kinetics 
process  as  the  flow  moves  downstream  from  n  to  n+1 .  The  sum  of  the  first  two 
terms  is  the  total  mass  fraction  of  category  particles  present  downstream  due 
to  kinetic  processes,  and  the  remaining  term  represents  the  contribution  to 
the  mass  fraction  from  diffusion.  Noting  Fig.  19,  the  diffusional  term  can 
be  split  into  the  contributions  from  each  of  the  three  upstream  points: 
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For  each  of  the  categories  present  at  a  given  grid  point,  we  may  express  the 
category  mass  fraction  as 


6(0,I)n+l  ,m 


+  aG£  +  AG3 


(22) 
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where  the  kinetics  contribution  is  evaluated  as  G(J,I)^  with  its  associated 
R(J,I)^  and  Knowing  the  category  sizes  associated  with  each  upstream 

contribution  the  Sauter  criteria  (or  other  appropriate  characterization 
criteria)  can  be  used  to  define  the  new  category  radius: 

R(J,I)n+l ,m  = 


G(J,I)k  +  AG1  +  AGg  +  AGg 
G(j  AG2  aG^ 

Rptlik  +  R(J,I)n,m+T  +  R(J,I)n,m  +  R(J,I)n,m-l 

Knowing  the  new  mass  fraction  and  size,  the  number  density  of  particles  is 
determined  from: 


(23) 


N(J,I)n+l  ,m 


pn+1 ,mG(J,I )n+l ,m 
4/3iip^R( J ,  I  )n+l  ,m 
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3.4  SUMMARY 

The  particle  tracking  and  boron  particle  combustion  subroutines  described 
in  this  section  represent  major  developments  with  respect  to  the  assembly  of  a 
modular  model  formulation  for  the  boron-slurry  fueled  ramjet  combustor. 

However,  considerable  further  work  still  remains  to  complete  the  modular  model 
development.  Some  of  this  work  involves  the  definition  of  improved  models 
for  the  kinetics  of  the  boron  combustion  process  and  thus  requires  additional 
experimental  data. 

As  was  noted  in  the  discussion  of  boron  particle  combustion  models  in 
this  section,  neither  of  the  models  presently  available  is  complete:  King's 
model  (Refs.  14,  15)  ignores  chemical  kinetic  effects  and  describes  only  the 
ignition  stage,  while  Edelman's  model  (Refs.  19,  20),  which  involves  kinetically 
controlled  combustion  phenomena,  ignores  the  oxide  layer  in  its  formulation. 

Of  the  two,  the  most  easily  extended  is  the  Edelman  model,  and  the  generaliza¬ 
tion  of  this  model  to  incorporate  the  effects  of  the  oxide  layer  on  the  combus¬ 
tion  process  is  presently  being  investigated. 

Further  development  with  respect  to  the  modeling  of  two-phase  flow 
phenomena  within  the  context  of  the  modular  formulation  involves  the 
incorporation  within  the  model  of  the  effects  of  thermal  nonequilibrium  and 
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modeling  of  the  simultaneous  oxidation  of  boron-containing  species  and  the 
hydrocarbon  species  in  the  gas  phase.  Work  in  these  areas  is  also  in  progress. 


4.  FORMULATION  OF  A  MATHEMATICAL  MODEL  OF  A  DUCTED  ROCKET  COMBUSTOR 


The  ducted  rocket  has  certain  features  that  are  closely  related  to 
ramjet  systems,  particularly  integral  rocket  ramjets  which  utilize  the 
sudden-expansion  ("dump")  combustor  concept  coupled  with  multiple  side  entry 
air  inlets.  The  ducted  rocket  is  sometimes  referred  to  as  a  gas  generator 
or  solid  propellant  rocket-fueled  ramjet.  Actually,  it  is  the  rocket  feature 
that  constitutes  the  major  differentiating  feature  between  ducted  rockets 
and  liquid  fueled  ramjets,  since  the  rocket  exhaust,  which  provides  the  fuel 
injection  into  the  ramburner,  is  choked  and  its  momentum  is  significant. 

Because  of  the  underexpansion  of  the  rocket  exhaust  local  supersonic  flow 
is  encountered  and  the  significant  level  of  the  rocket  exhaust  momentum 
relative  to  the  airflow  momentum  can  strongly  affect  the  overall  ramburner 
flowfield.  These  factors  are  not  present  in  a  liquid  fueled  ramjet.  Signifi¬ 
cant  differences  also  exist  between  the  ducted  rocket  fuel,  which  is  made  up 
of  the  products  of  combustion  of  a  fuel-rich  solid  propellant,  and  the  liquid 
(or  slurry)  fuels  used  in  ramjets.  Further,  while  practical  sudden  expansion 
ramjet  combustors  can  be  axi symmetri c ,  ducted  rocket  configurations,  as  a 
rule,  are  three-dimensional.  Thus  the  analytical  simplifications  offered 
by  axisymmetry  are  not  available  in  modeling  the  ducted  rocket  flowfield. 

Therefore,  although  there  is  a  logical  connection  between  the  combustion 
processes  in  ramjets  and  in  ducted  rockets,  it  is  clear  that  there  are 
essential  differences  that  must  be  addressed  in  order  to  establish  a  better 
quantitative  understanding  of  the  mechanisms  and  parameters  involved  in  the 
ducted  rocket  combustion  process.  Despite  the  evident  difficulties,  develop¬ 
ment  problems  that  have  been  encountered  clearly  show  the  need  for  research 
into  ducted  rocket  combustor  aerodynamics:  flow  patterns,  mixing  rates, 
and  chemical  kinetic  rates  in  ducted  rocket  environments  need  to  be  defined, 
and  analyses  tailored  to  the  specific  characteristics  of  the  ducted  rocket 
combustion  chamber  need  to  be  developed. 

The  complexity  inherent  in  the  ducted  rocket  combustion  chamber  argues 
for  the  use  of  a  modular  formulation  such  as  has  been  successfully  developed 
for  axisymmetric  liquid  fueled  sudden  expansion  ramjet  combustor  configurations. 
However,  while  a  data  base  which  provides  some  description  of  the  combustor 
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aerodynamics  exists  for  the  dump  combustor,  no  such  data  base  exists  in  the 
case  of  the  ducted  rocket.  Since  the  modular  approach  depends  on  an  a  priori 
description  of  the  size  and  location  of  major  flowfield  elements,  particularly 
regions  of  recirculating  flow,  an  initial  data  base  is  critical  to  the 
development  of  a  detailed  modular  model.  In  the  absence  of  experimental 
data,  some  information  on  the  aerodynamics  involved  in  the  ducted  rocket 
combustion  chamber  flowfield  can  be  obtained  through  the  use  of  a  three- 
dimensional  Navier-Stokes  solution  procedure.  However,  this  type  of  modeling 
is  complex  and  not  readily  amenable  to  the  incorporation  of  the  detailed 
chemical  kinetics  which  is  also  required  to  gain  an  overall  analytical  descrip¬ 
tion  of  ducted  rocket  flowfield  details.  It  is,  on  the  other  hand,  possible 
to  use  the  Navier-Stokes  solution  procedure  either  for  nonreacting  flow  or, 
using  a  simple  representation  of  the  primary  heat  release  reactions,  for 
reacting  flow,  to  provide  the  detailed  information  necessary  to  construct  an 
appropriate  modular  model.  This  is  the  approach  adopted  in  the  present  work. 

In  order  to  discuss  the  requirements  on  a  three-dimensional  Navier-Stokes 
solution  procedure  for  the  ducted  rocket  it  is  useful  to  consider  a  specific 
configuration.  To  this  end,  the  configuration  which  is  to  be  utilized  in  an 
experimental  program  at  the  Aeropropulsion  Laboratory  (AFWAL)  is  of  interest: 
a  sketch  of  this  configuration  is  shown  in  Fig.  ZO.  Basic  modeling  requirements 
for  the  main  combustor  include  the  use  of  the  three-dimensional  form  of  the 
Navier-Stokes  equations  in  cylindrical  coordinates,  along  with  a  suitable 
turbulence  model  and  the  capability  to  provide  the  particle  tracking  required 
to  compute  two-phase  flow  transport  phenomena  as  outlined  in  the  previous 
section.  For  this  specific  configuration  a  half-domain  computation  is 
sufficient  but  for  generality  the  capability  to  compute  the  full  flowfield 
domain  may  be  necessary. 

Two  significant  complications  are  inherent  with  this  flowfield.  First, 
the  air  inlet  to  the  combustor  is  provided  through  rectangular  ducts,  so  that 
the  inlet  airflow  is  quite  nonuniform,  and  second,  the  flowfield  in  the  region 
of  the  gas  generator  nozzles  must  also  be  computed.  In  the  latter  case,  the 
scale  of  the  flow  is  significantly  smaller  than  the  scale  which  characterizes 
the  overall  flowfield.  To  adequately  handle  both  of  these  complications, 
variable  spacing  capability  is  required  using  both  nonuniform  mesh  spacing  and 
(possibly)  coordinate  transformations. 
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FIGURE  20.  Ducted  rocket  combustor  configuration, 
AFWAL/APL. 
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The  flowfield  within  the  inlet  ducts  involves,  since  recirculations  are 
not  expected  to  be  present,  a  three-dimensional  parabolic  or  quasiparabolic 
flow.  In  general,  the  flow  in  this  duct  would  be  coupled  to  the  conditions 
within  the  combustor,  but  at  least  as  a  first  approximation  an  uncoupled 
analysis  can  be  used.  One  aspect  of  this  flowfield  that  requires  special 
attention  is  the  matching  of  the  inlet  duct  grid  with  that  utilized  for  the 
main  combustor,  as  shown  schematically  in  Fig.  21.  Matching  of  grids  is 
also  a  requirement  with  respect  to  the  injector  nozzle  face,  as  also  shown  on 
Fig.  21.  The  most  suitable  approach  to  this  aspect  of  the  problem,  given  the 
different  scales  which  characterize  the  injector  nozzle  and  the  combustor  as 
a  whole  is  the  use  of  a  control  volume  analysis  for  the  computational  region 
surrounding  the  injector  nozzles. 

Development  of  an  aerodynamic  model  for  the  combustor  shown  schematically 
in  Fig.  20  has  been  based,  initially,  on  the  use  of  the  three-dimensional 
elliptic  formulation  developed  at  Imperial  College,  and  described,  for  example, 
by  Patankar  (Ref.  21).  This  code  is  a  primitive-variables  approach,  which  is 
to  say  that  it  solves  the  equations  for  the  mean  velocity  components  directly, 
rather  than  using  transformed  variables  such  as  the  stream  function  and 
vorticity.  The  solution  is  obtained  iteratively  using  an  underrelaxation 
technique,  with  a  semi-implicit  pressure-correction  procedure  and  a  velocity 
update  scheme.  For  the  ducted  rocket  application  the  code  requires  considerable 
modification,  and  during  this  modification  effort  the  structure  of  the 
technique  is  being  revised  to  allow  modifications  to  various  aspects  of  the 
model  (i.e.,  turbulence  models,  two-phase  flow  formulations,  and  simple 
chemical  kinetics  approaches)  to  be  more  readily  incorporated.  The  basic  code 
is,  however,  a  steady  state  formulation;  an  unsteady-state  approach  may 
eventually  be  required  because  of  the  oscillatory  nature  of  the  flow  that  has 
been  observed  in  some  ducted  rocket  combustor  experiments. 
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FIGURE  21.  Computational  grid  requirement  for  ducted  rocket 
combustor  computation. 


5.  OVERALL  SUMMARY  AND  FUTURE  WORK 


Research  effort  for  this  period  has  focused  on  three  major  areas  of 
investigation:  further  examination  of  the  use  of  the  modular  ramjet  combustor 
model  as  a  technique  for  the  analysis  and  interpretation  of  combustor  test 
data,  the  development  and  implementation  of  a  modular  model  for  the  boron 
slurry  fueled  ramjet,  and  the  development  of  analytical  models  for  ducted 
rocket  combustors.  Further  examination  of  the  liquid-fueled  ramjet  modular 
model  has  shown  that  the  sti rred-reactor  representation  of  the  flame-stabilizing 
recirculation  zones  in  a  dump  combustor  provides  a  useful  and  powerful  tool 
for  the  examination  of  dump  combustor  flame  stability  characteristics. 

This  work  has  also  demonstrated  the  utility  of  the  model  in  terms  of  the 
analysis  and  interpretation  of  ramjet  test  data.  One  area  of  required  future 
work  has  been  shown  to  involve  the  engineering  modeling  of  chemical  kinetics 
phenomena  under  fuel -rich  conditions:  in  some  ramjet  designs  and  under  some 
operating  conditions  local  regions  within  the  combustor,  and  particularly 
the  recirculation  regions,  can  be  considerably  more  fuel-rich  than  would  be 
expected  based  on  the  overall  fuel -air  ratio  for  the  operating  condition. 
Engineering  chemical  kinetics  models  for  fuel -rich  conditions  are  being 
developed  by  SAI  under  other  research  programs,  and  future  work  in  this  area 
of  the  current  AFOSR  program  will  involve  the  incorporation  of  these  models 
into  the  modular  formulation  of  the  liquid-fueled  sudden-expansion  ramjet. 

Initial  work  relevant  to  the  boron-slurry  fueled  ramjet  has  centered 
on  the  analysis  and  implementation  of  existing  models  for  boron  particle 
combustion  and  the  development  and  implementation  of  a  particle  tracking 
technique  for  use  with  the  modular  model  formulation.  Existing  models  for 
boron  particle  ignition  and  combustion  do  not  encompass  the  full  range  of 
combustion  phenomena  relevant  to  boron  combustion  necessary  to  provide  a 
comprehensive  model.  Future  work  in  this  area  involves  the  extension  of 
existing  approaches  to  consider  both  diffusion-  and  chemical  kinetically- 
controlled  combustion  phenomena,  including  the  effects  of  an  oxide  surface 
layer  on  the  particle.  The  particle  tracking  methodology  developed  for  use  in 
the  modular  model  provides  the  only  realistic  method  of  computing  the  combus¬ 
tion  behavior  of  particle  clouds,  but  this  approach  needs  to  be  extended  to 
include  thermal  nonequilibrium  between  the  particles  and  gas  phase  and  the 
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simultaneous  combustion  of  boron  and  gas-phase  hydrocarbon  species.  Future 
work  will  involve  these  extensions. 

An  aerodynamic  model  for  the  ducted  rocket  combustor  configuration  is 
under  development,  based  initially  on  the  three-dimensional  elliptic  formula¬ 
tion  developed  at  Imperial  College.  Computations  of  the  flowfield  within  a 
representative  ducted  rocket  configuration  will  be  carried  out  during  the 
next  year's  work.  The  results  of  these  computations  will  be  used  to  assess 
the  utility  of  the  approach  (including  the  possibility  of  there  being  a 
requirement  for  the  development  of  a  time-dependent  model)  and  to  provide  a 
baseline  of  flowfield  detail  required  to  initiate  the  development  of  a  detailed 
modular  model  of  this  flowfield. 
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